% Compare average before and after 
% passing the 2nd order gas phase
for i = 1:20
    fprintf(1,'i = %d, --------------\n',i);
    w1 = 0.4+0.04*(rand(1000,1)-1);
    % Path#1: average before (InvReaction + ReactionModel)
    % ReactionModel(InvReactionModel(mean(w1))
    w1_mean = mean(w1);
    x1_mean = InvReactionModel(w1_mean,133.32,850);
    w1_mean_recover(i) = ReactionModel(x1_mean,133.32,850);
    w1_mean_recover(i)
    % Path#2: average after (InvReaction+ReactionModel)
    % mean(ReactionModel(InvReactionModel(w1)))

    x1 = InvReactionModel(w1,133.32,850);
    w1_recover = ReactionModel(x1,133.32,850);
    w1_recover_mean(i) = mean(w1_recover);
    w1_recover_mean(i)
end
figure(1);
hist(w1_mean_recover-w1_recover_mean);
title('w1.mean.recover-w1.recover.mean');


